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Abstract 

Lepton asymmetry generated in the early universe by neutrino oscillations 
into sterile partners is calculated. Kinetic equations are analytically reduced 
to a simple form that permits an easy numerical treatment. Asymptotic values 
of the asymmetry are at the level of 0.2-0.3 and are reasonably close to those 
obtained by other groups, though the approach to asymptotics in some cases is 



^ . noticeably slower. No chaoticity is observed. 

2 : 1 Introduction. 



Neutrino oscillations in the early universe might possess a very interesting property, 



Qh! if active neutrinos {ua = I'e, i^fj., ^t) are mixed with sterile ones (z/^)- Refraction index 
'-^ ' of neutrinos in the cosmic plasma depends upon the cosmological charge asymmetry 

of the plasma, rj, which is normally quite small, rj = 10~^ — 10^^'^. Because of this 
dependence the transformation of Ua into Us might be slightly more favorable than the 
transformation of the corresponding antineutrinos (or vice versa depending upon the 
sign of the initial asymmetry). The feedback effect is positive and leads to a further 
increase of asymmetry making, say, Ua Vs transformation more and more efficient 
in comparison with z/^ — > z/^. Of course the total leptonic charge of active plus sterile 
neutrinos is conserved but refraction index depends only on charge asymmetry in the 
sector of active neutrinos, and the lepton asymmetry in active sector could strongly 
increase. 



lAlso: ITEP, Bol. Cheremushkinskaya 25, Moscow 113259, Russia. 
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The effect of asymmetry generation takes place only for a sufficiently small mix- 
ing; i.e. for 5m? ~ 1 eV^ it is roughly bounded by (sin 26')^ < (~ 10~^), where 
9 is the vacuum mixing angle. For large mixings active and sterile states of both 
neutrinos and antineutrinos would quickly reach thermally equilibrium values, they 
would become equally populated, and this prevents from generation of a large asym- 
metry. For a positive mass difference between Vs and Va {5m? > 0) the asymmetry 
would remain small also in the case of small mixing. However for < the reso- 
nance MSW-transition Q might take place in the primeval plasma and this effect 
could compensate a smallness of vacuum mixing and produce a considerable lepton 
asymmetry in the sector of active neutrinos. 

The instability with respect to generation of lepton asymmetry by neutrino oscil- 
lations was noticed in ref. but it was concluded there, on the basis of simplified 
arguments, that the rise of asymmetry was terminated when it was still quite small. 
This conclusion was reconsidered in ref. |Q (see also refs. where it was argued 

that a very large asymmetry, close to 1, may be generated by the oscillations. This 
result was questioned in ref. [^j where it was claimed that back reaction effects are 
very strong, they significantly slowed down the rise of asymmetry, and the latter 
could reach only the value around 10"^ — 10~^ for reasonably small mass differences. 
However certain drawbacks of the approach of the paper were indicated in ref. |10| 



(see also |[TT| , |T2|) and though not all of them were relevant, it would be enough to 
have one weak point to destroy a conclusion. 

Since the result of a large asymmetry generation is very interesting and impor- 
tant (in particular, for the big bang nucleosynthesis, BBN) it is worthwhile to make 
an independent calculation of the effect. In this work kinetic equations governing 
neutrino oscillations in the early universe are analytically transformed to a rather 
simple form which after well controlled approximate simplifications permits an easy 
numerical solution. 

The results of this work are reasonably close to those found in refs. 0, |6|, [^, p!0| 
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for the case of {i^fi^T — i^s)-oscillations and approximately 2-3 times smaller than the 
asymmetry found in ref. |]l3l for (z/g— z^s)-oscillations in the temperature range essential 
for BBN. 

Another interesting issue related to the asymmetry generation is its possible 
chaoticity so that the sign of the final large value of the lepton asymmetry is a rapidly 
oscillating function of the oscillation parameters or of a small variation of the initial 



value of the cosmological charge asymmetry |jT^-||2T|. According to the calculations 
of the present paper, no chaoticity exists at least in the parameter range considered 
below, for which a large lepton asymmetry is generated. 

2 Kinetic equations. 

We assume that a non-negligible mixing exists only between one active and one sterile 
neutrino, so that the neutrino state is described by 2 x 2 - density matrix p. Its 
evolution is governed by the usual equation: 

tp = [n,p] (1) 

where Ti is the neutrino Hamiltonian. It contains the free part: which is diagonal in 
the mass basis: 



Hfree = diag 



p'^ + ml, Jp^ + m 



(2) 



and the part that describes interaction with medium which looks simpler in the flavor 
basis. The interaction Hamiltonian contains first order terms given by refraction 
index (or effective potential of neutrinos). It was originally calculated in ref. and 
consists of two terms (see also discussion in papers p^, ||: 

V^ff = ±Ci7]^''^GfT^ + C^^^^ , (3) 



a 
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where E is the neutrino energy, T is the temperature of the plasma, Gp = 1.166 ■ 10 
GeV~^ is the Fermi coupling constant, a = 1/137 is the fine structure constant, and 



the signs "±" refer to anti-neutrinos and neutrinos respectively. According to ref. p2 
the coefficients Cj are: Ci ~ 0.95, C| ~ 0.61 and C!^'^ ~ 0.17. These values are true 
in the case of thermal equilibrium, otherwise these coefficients are some integrals over 
the distribution functions. The contributions to 77*^"^ from different particle species 
are the following: 

ri^^'^ = 2r]^^ + rj^^ + rj^^ + rj^ - r/„/2 (for u^) , (4) 
r/(^) = 2ri,^ + ri,^ + r],^ - r]n/2 (for u^) , (5) 

and rf^^^ for is obtained from eq. (|^) by the interchange ji ■>r^ r. The individual 
charge asymmetries, rjx-, are defined as the ratio of the difference between particle- 
antiparticle number densities to the number density of photons: 

^x = {Nx-Nx)/N, (6) 

The interaction Hamiltonian contains also the so called second order terms that 
describe the loss of coherence due to elastic or inelastic neutrino scattering and anni- 
hilation as well as neutrino production by reactions in the primeval plasma. For the 
description of the loss of coherence it is enough to take into account only imaginary 
part of the second order Hamiltonian. Discussion and calculations can be found in 



the papers |2^-[0. The exact form of the second order terms has quite complicated 
matrix structure, it is non-linear in p, and contains multi-dimensional integrals over 
phase space (which can be reduced down to two dimensions). Their explicit expres- 
sions can be found e.g. in refs. 29 1. However it is very difficult, to solve kinetic 



equations for density matrix with the exact expressions for the second order terms, 
especially in the resonance case. (An approximate solution in non-resonance case is 
found in ref. |3^). Hence one usually mimic the exact expression by the linear "poor 



man" substitution: 

^p = ...(V2){r, p-pej... (7) 

where curly brackets mean anti-commutator and the matrix F effectively describes 
neutrino scattering and annihilation. It has the only non- vanishing entry at aa-corner 



in the flavor basis, which is taken as 

r,, = = c^GlT^E (8) 

where is a constant. According to reference Cp = 1.27 and C^^'^^ = 0.92, while 
according to 0, where slightly more accurate calculations were done, = 1.13 and 
_ Q YQ. This difference and even the absolute value of Taa are not essential for 
the magnitude of rising asymmetry because F remains in some sense small and can be 
neglected (see below). Even the use of the approximate expression for the coherence 
breaking terms (0) in contrast to the exact collision integral which changes the results 



in the case of non-resonance oscillations by an order of magnitude |Q, have a very 
weak impact on the value of the generated lepton asymmetry in the resonance case. 

The matrix peq is equal to feq{lj) I where / is the unit matrix and and feq{E/T, fj,) 
is the equilibrium Fermi distribution function: 

^""^ ^ exp [{E - /i) /T] + 1 • 
where /i is neutrino chemical potential. In what follows we will take /i = 0, though in 



refs. @]-[§,||T0| it was taken equal to the running value determined by the generated 
lepton asymmetry. However the final result for the asymmetry is not sensitive to this 
assumption (see discussion in section ^. 

In the Friedman-Robertson- Walker universe the matrix elements of p satisfy the 
following equations: 

i{dt- Hpdp)paa = Fo{psa- pas)/'^-i^o{paa- feq) , (10) 

i{dt- Hpdp)pss = -Fo{psa- pas)/^ , (11) 

i{dt - Hpdp)pas = WoPas + Fo{pss - Paa)/'2 - «FoPas/2 , (12) 
i{dt - Hpdp)psa = -WoPsa - Eo{pss - Paa)/2 - i^o/2psa , (13) 

where a and s mean "active" and "sterile" respectively, 

Fo = W sin 2^/2E, (14) 
Wo = 5m^ cos 20 /2E + l/^p (15) 
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H = ySnptot/'^Mpi is the Hubble parameter, Mpi is the Planck mass, and p is the 
neutrino momentum. The neutrinos are assumed to be very light, so that p = E. The 
total cosmological energy density is taken as ptot = 10.757r^T^/30. This corresponds 
to photons, e^'^-pairs, and three types of neutrinos in thermal equilibrium with equal 
temperatures T. This approximation is quite good when T is larger than a few 
hundred keV. At smaller temperatures e"'"e~-annihilation heats up electromagnetic 
component of the plasma and the temperatures of photons and neutrinos become 
different. 

We have prescribed sub-index "0" to the coefficient functions in eqs. (|10|-|13]) to 
distinguish them from the same ones divided by Hx (see eqs. ( p7| - |3TD below). 



The anti-neutrino density matrix, p, satisfies the similar set of equations with 
the opposite sign of the antisymmetric term in V^jj- and with a slight difference in 
damping factors Tq which is proportional to the lepton asymmetry. 

It is convenient to introduce new variables: 

X = moR{t) and y = pR(t) , (16) 

where R(t) is the cosmological scale factor so that H = R/R and mg is an arbitrary 
mass (just normalization), taken as tuq = 1 MeV. One may approximately assume 
that T = —HT, and correspondingly R = 1/T. In terms of these variables the 
differential operator {dt — Hpdp) transforms to Rxd^- We will normalize the density 
matrix elements to the equilibrium distribution with zero chemical potential: 

Paa = feq{y) [1 + ^{x, ?/)], Pss = feqiv) [1 + s{x, ?/)] , (17) 

pas = psa = feq{y)[h{x,y) + il{x,y)] , (18) 

Other authors find it convenient to express the density matrix in terms of Pauli 
matrices and the polarization vector, P = {P^, Py, Pz), in such a way that: 



l+P-a , (19) 
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The relation between different functions are PqPx = 2feqh, P^Py = —2feql, PqPz = 
feq{a. — s) and Pq = feq{2 + a + s). In particular, Pz = 1 means that all the neutrinos 
are active, z/g, t'^, or Ur- 



Let us now divide both sides of equations (|T0|-|T3D by Hx and denote the cor- 
responding coefficient functions ( p!^ - |T5|) with sub-index "1", e.g. Wi = Wq/Hx = 
Ui ± ViZ, etc and thus we get: 

„2 



z 



1.12 ■ 109cos2^(5m2— + ^ 
y ^ x^ 

30 

~2 ' 



X 

10^0 
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dy 



y^feqiy) (a - a] 



where is expressed in eV^ (here and below), C*! 



26 and C^'" 



(20) 
(21) 
(22) 

7 and rjo 



is the charge asymmetry of all particles except for z/^ defined in accordance with 
eqs. (§J^). The normalization of the charge asymmetry term ( ^2]) is rather unusual 
and to understand the numerical values of the coefficients one should keep in mind 
the following. The coefficient C2 in eq. (|^) is found for the standard normalization 
of charge asymmetry with respect to the present day photon number density which 
differs from that in the early universe by the well known factor 11/4, related to 
the increase of photon number by e+e~-annihilation. On the other hand, lepton 
asymmetry, L^^, induced by neutrino oscillations, which is calculated in most of 
papers is normalized to the number density of photons that are in thermal equilibrium 
with neutrinos, so the factor 11/4 is absent. The photon number density is equal to 
= 2({3)T^/n'^ (C(3) ~ 1-2). The charge asymmetry of neutrinos is 

Vu = J dyy'^ipaa - Paa) (23) 

so that T]^^ = 4Ly^/ll. The quantity Z introduced in eq. 
factor 2 ■ 10-i°7rVC(3): 



differs from L^^ by the 



L = 16.45 ■ 10"^° Z 



(24) 
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(25) 



The factor 10 is chosen so that initially Z = 0(1). Noting that the charge asym- 
metry of neutrinos under study enters expressions (§,|^) with coefficient 2 one obtains 
the coefficient 1/11.96 ^ 1/12 in eq. (|22|). 

Now, following ref . , we will introduce one more new variable q = ^^x^ in such a 
way that Ui, eq. (pO]), vanishes at q = y independently of the values of the oscillation 
parameters: 

2 

Ui = 1.12-10^ cos29 \5m'^\ — 

y 

(it is assumed that < 0). The coefficients are 

^e = 6.63-10^(|5m2|cos2^)^^\ = 1.257 ■ 10^ (l^m^j cos2^)^^^ (26) 

In what follows we assume that sin 20^1 and cos 26* ^ 1 . 

Written in terms of the variable q the system of basic kinetic equations takes a 
very simple form 

s' = -{Kjy)sm2ei (27) 

a' = {Ka/y){sm29l-2-fa) (28) 

h' = {Kjy){Wl-^h) (29) 

I' = {Ka/y)[sm2e{s - a)/2-Wh~^l] (30) 

where Ka = 1.12-109 cos20 |(5m2|/3^„, so = 5.63 ■ 10^(|(5m2| cos 2^)^/2 and iT^,^ = 
2.97 ■ 10^(|5m2| cos 2^)1/2. In what follows we use the limit Ka ^ 1. This permits 
to make accurate analytic calculations. A large magnitude of this coefficient reflects 
a large frequency of neutrino oscillations with respect to other essential time scales. 
Its large value makes numerical solution very difficult but allows to make accurate 
analytical estimates. 

The coefficient functions in these equations have the form: 

W = U±yVZ, U = y^q-'^-l, 

V = Kq-'l\ J = eay'q-' (31) 



where the signs " or in W refer to neutrinos and antineutrinos respectively; 
be = 3.3 ■ 10-3(|5m2| cos2^)-i/^ 6^,^ = 7.8 ■ IQ-^d^m^l cos2^)-i/3^ and ea are small 
coefficients, eg ~ 7.4 • 10^'^ and e^^r ~ 5.2 ■ 10~^. Their exact numerical value is not 
important. It is noteworthy that the charge asymmetric term in W comes with a 
very large coefficient if expressed in terms of L, VZ ~ 10''g~^/^L, while in all other 
possible places L (or chemical potential) enters with the coefficient of order 1. 
It follows from eqs. (I^-IOD that 



dg (a' + s^ + 2h^ + 2f) = -4-f{K/y) (a^ + + f) (32) 
so that the quantity in the r.h.s. may only decrease. 

3 Solution of kinetic equations. 

One can solve analytically the last two kinetic equations (|^,^) with respect to h 
and / in terms of a and s: 

liq, y) = -{K sin 29 /2y) f dq, [a(gi) - s(gi)] e"^^ cos A$, (33) 
h{q, y) = -{K sin 26 /2y) f dq, [a{q,) - s{q,)] e"^^ sin A$, (34) 

where qm is the initial "moment" q from which the system started to evolve, AF = 
r(g, y) - T{qi,y), A<l> = <^{q, y) - <i>{qi,y), and 

d,r = K^/y, d,^ = KW/y. (35) 

We rewrite the ffist two equations (P7| , P5| ) in terms of a = a + s and 5 = a — s: 

a' = ~{K^/y){a + 5) (36) 
5' = {2Ksm2e/y)l-{K'y/y){a + S) (37) 

The ffist of these equations can be solved for a: 

cr{q, y) = aUy) e-r('''^)+^"(^) - - T dq^e'''^^{q^,y)5{q,, y) (38) 



y Jq. 
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The first term in tliis expression, proportional, to the initial value am is exponentially 
quickly "forgotten" and we obtain the following equation that contains only 5 (and 
another unknown function, integrated charge asymmetry Z{q) that is hidden in the 
phase factor A$): 



y 

'Ksin 29^ 

. y . 



6 + 



dqie ^^l{qi,y)5{quy) 



/ rfgi(5(gi,?/)e"^^cosA$ 



(39) 



(40) 



Up to this point this is an exact equation (with the omitted initial value of a which 
contribution is exponentially small). There is also an uncertainty related to the choice 
of the form of peq in eq. (0) either with zero or non-zero chemical potential, see eq. (P). 
We have chosen here /i = and, as is argued in sec. |^, the choice fi does not lead 
to a noticeably different results. This ambiguity could be rigorously resolved if one 
uses exact collision integrals instead of eq. (0). This will be discussed elsewhere. 

Let us take now similar equation for antineutrinos and consider the sum and dif- 
ference of these two equations for charge symmetric and antisymmetric combinations 
of the elements of density matrix, = 6 + 6 and A = 6 — 6. The equations have the 
following form: 



y2 

fKsin 26^ 



dqi e 



-Ar 



(41) 



y2 

fKsm 26^ 



/ dqi e"^^7i Si 



dqi e 



-Ar 



(42) 



where sub-1 means that the function is taken at gi, e.g. 71 = 7(^1, y), etc; c = cos A$, 



and c = cos A$. Using expressions (|3l| , |35|) we find: 



c — c 
2 

c + c 



sm 



cos 



\ y qqi^ 

K{q-qi)(--+ y ^ 



y 



10 



sm 



cos 



K r dq2V{q2)Z{q2) 
K r dq2V{q2)Z{q2) 



(43) 
(44) 



Here V{q) is given by expression (|3l|) and does not depend on y. 

At this stage we will do some approximations to solve the system ( ^^21) . First, 



let us consider the terms proportional to 7. They are definitely not important at 
large q (or low temperature). Let us estimate how essential are they at low g (g ~ 1). 
Integrating by parts the first term in the r.h.s. of eq. (0), using expression (|35D and 
neglecting exponentially small contribution of the initial value, we find: 

A — / dqie 71^1 = / dqic —— (45) 

y y Jq,n y Jq,n «gi 

The remaining integral can be easily evaluated in the limit of large K'-f/y = Key/q^. 
Indeed, AF = Kty{q — qi)/qqi and for q <1 the coefficient in front of the exponential 
(g — gi) is larger than 400 for v^. and 300 for and v^-. So the integral strongly sits 
on the upper limit and, together with the coefficient in front, it gives just 5'{q). 
Thus, when the 7-terms are large they simply double the coefficient of A' in eq. (|^): 
A' 2A'. A possible loophole in this argument is a very strong variation of the 
integrand, much stronger than that given by exp(AF). However one can check from 
the solution found below that this is not the case. 



Thus the role of 7 terms in eq. (41) is rather mild, they could only change the 
coefficient in front of A' from 1 to 2, and become negligible for large g where the 
bulk of asymmetry is generated. So let us neglect these terms in the equation. This 
simplification does not have a strong impact on the solution. 

Let us make one more approximate assumption, namely let us neglect the second 
term, proportional to Ai, in the last integral of the r.h.s. of eq. (0). Initially 
S = 2 and A = 10~^ — 10~^° and the neglect of A in comparison with S is a 
good approximation, at least at initial stage. We will check the validity of this 
assumption after we find the solution. And last, we assume that S changes very 
slowly S Sj„ = 2. Justification for the latter is a smallness of the mixing angle, 
sin 2^ ~ 10^"^. In the limit of zero mixing, the solution of eq. (42) is E = const. We 



will relax both these assumptions in sec. ^. 

As a last step we need to find a relation between A = a — s — a + s and charge 
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asymmetry Z . To this end one may use the conservation of the total leptonic charge: 

/•oo 

/ dyy^ feq{y) {a + s - a - s) = const (46) 
Jo 

Using this conservation law we find: 
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10 



d_ 

dq 



,dZ 
dq 



(47) 



Keeping all these assumptions in mind we can integrate both sides of eq. ( ^Tj ) with 
dyy'^feqiy) and obtain a closed ordinary differential equation for the asymmetry Z[q), 
valid in the limit of large K. Integration over y gives: 



in" 



dyfeqiy) 



dqi exp 



sm 



\y QQiJ 



sm 



bK I dgs^Tg- 
<?i g2 



(4J 



where b is defined in eq. (^Tj) and C = K{q — qi). Integration over y here can be done 
explicitly and the result is expressed through a real part of a sum of certain Bessel 
functions of complex arguments. To do that one has to expand 



/e, = E(-ir^'exp(-ny) 



(49) 



and integrate each term analytically [^]. It can be seen from the result (it is more 
or less evident anyhow) that the integral over qi is saturated in the region ^ ~ 1. So 
that we can take qqi ~ and 



K rdq2Z{q,)q^'/'^CZ{q)/q'^' 

Jqi 



(50) 



The correction in this expression is of the order of Z'{q)/K. It can be checked, using 
the solution obtained below, that the correction terms are indeed small. 



Keeping this in mind we can take the integral over ( in the r.h.s. of eq. (HSI) . To 
ensure convergence we proceed as follows. First, we expand feq{y) in accordance with 
eq. (|49|). Each term of the series is non-singular in the complex y-plane. After that 
we can rotate the contour in the complex y-plane to imaginary axis, clockwise for the 
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exp[i^(l/?/ — ?//g^)]-part of sin[C(l/y — ylq^)], {y — > —iy), and counter-clockwise for 
the complex conjugate part {y iy). Both terms give exp[—({l/y + y/q^)], so the 
integral over C, is exponentially converging for each term in the series (|i9| ) and the 
resulting series is convergent as well: 

10i°fs:(sin2^)2 .CO 



^^^^ - 8^^^ 4^^"^^ Jo 47 



n 

1 



+e 



-my 



1 

. y 9"' ^ y ^ 



(51) 



where = = bZ(q)q~^^^. 

Since e ~ 10~^ is a small number it may be neglected and changing the integration 
variable, y = qt, we come to the expression: 

where x(Q') = q'Piq) = bZ{q)q^^^^. Both the integral over t and summation over n 
can be done explicitly and we finally obtain: 

-Z' = ^°'°^/g7^^^' ^^f^^ [t2Uqh)-tMqh)] (53) 

where we introduced the coefficient n, such that k = 1 for large q and k, = 2 for 
g ~ 1. It reflects the role of decoherence terms, proportional to 7, see discussion 
after eq. (^Sf). The quantities ti 2 are the poles of the denominator in eq. ( ^2]) in the 
complex upper half-plane of t (resonances): 

^1,2 = ^ (54) 

It is an ordinary differential equation that can be easily integrated numerically. It 
quite accurately describes evolution of the lepton asymmetry in the limit when back 
reaction may be neglected: we assumed above that S = 2 and A <^ S. 
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Before doing numerical integration let us consider two limiting cases of q close to 
initial value when asymmetry is very small and the case of large q. When q is not 
too large, g ~ 1, the r.h.s. of eq. (|^) can be expanded in powers of x and we obtain 
a very simple differential equation that can be integrated analytically: 

^, 10i0ir(sin2e)2 2^ , ^ , ^r 1 , , f ^^^ /..n 

= — -q feMx{q) [-1 + 9 (1 - feM))] (55) 

(we took here k = 2). One sees that for q < qmm = 1-278 the asymmetry exponentially 
decreases and reaches the minimum value 



V 



exp 



(56) 



647r JO V l + exp(— g)^ 

The integral in the expression above is equal to 0.07539 and e.g. for (z/g — Vs)- 
oscillations the initial asymmetry drops by 3 orders of magnitude in the minimum. 
The drop would be significantly stronger even with a mild increase of mixing angle 
or mass difference. The temperature, when the minimum is reached (corresponding 
to qmin = 1.278) is 

T^,„ = 17.3 {Sm'f^' MeV, T^f, = 23.25 {Sm'f^' MeV (57) 

These results rather well agree with ref. []T3[ for Uf., while agreement for and 



case (see e.g. the papers ^ ^ |lOl) is not so good. 

For q > qmin the asymmetry started to rise exponentially and this regime lasted 
till X becomes larger than one and the asymmetry reaches the magnitude Z ~ 10'^ or 
L ~ 10^^. After that the asymmetry started to rise as a power of q. For large q and 
X the term containing ^2 dominates the r.h.s. of eq ( |55D and now it takes the form: 



Z^Z'.I^^^J^^^'WZ) ,58) 

where V is given by eq. (0). Assuming that VZ is a slowly varying function of q we 
can integrate this equation and obtain: 

Z{q) ^ 1.6 ■ lO^g^^/^ or L{q) ^ 2.5 ■ 10~^g"/^ (59) 
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The concrete values of numerical coefficients above are taken for (z/g — z^s)-oscillations 
with 5m^ = — leV^ and sin 2^ = 10"'*. This result is in a good agreement with the 
numerical solution of eq. ( |53D and the functional dependence, L ~ gi^/^ ~ j^-ii/s^ 



agrees with that found in ref. |]TT| and slightly disagrees with the results of refs. 0- 
0],[|1O|, |l^] where the law L ~ g*/^ ~ was advocated. 

Numerical solution of eq. (^) is straightforward. It well agrees with the simple 
analysis presented above. In the power law regime, where the bulk of asymmetry is 
generated, it is accurately approximated by the found above law (^): 

Le = 2.5- lO-^Ceg^^/^ (60) 

For Uf, — Vg mixing with (sin 20)^ = 10~® and 5m^ = —1 the correction coefficient 
Ce is 0.96, 0.98, 1, 1.01, and 0.997 for q = 6630, 1000, 100, 10, and 5 respectively. 
The results of numerical solution well agree with those of ref. in the temperature 
range from 10 down to 1 MeV for z/g — case. 

For the (z/^ — z/<j)-mixing with = —10 and (sin 2^)^ = 10~^ the solution can 
be approximated as 

L^^1.2-10"^C^g"/^ (61) 



with the correction coefficient = 0.84, 0.9, 0.98, 1.02, 1.05, 1 for g = 4 ■ 10^ 10^ 
10^, 10^, 10, 5 respectively. These results reasonably well agree with the calculations 
of ref. 1^ in the temperature range from 25 down to 2 MeV. At smaller temperatures 
this power law generation of asymmetry must stop and it abruptly does, according to 
the results of the quoted papers, but the solutions of eq. ([53| ) continue rising because 
back reaction effects are neglected there. We will consider these effects in the next 
section. 

4 Back reaction. 

The solution obtained above should be close to the exact one if A{q,y) <^ T,{q,y) 
and S ^ 2 = const, see eq. (^Tf). Since now we know the function Z{q) we can find 
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A(g, y) and check when this assumptions are correct. We will consider the region of 
sufficiently large q when the second pole (resonance) ti in eq. ( ^3]) is not important. Its 
contribution is suppressed as exp(— g^FZ) and it may be neglected already at g > 5 
(we assume for definiteness that initial value of the asymmetry is positive, otherwise 
the role of the two poles would interchange). In terms of oscillating coefficients cos A$ 
or cos A$, entering eq. (l4ll) , it means that only one of them is essential. It has a saddle 
point where the oscillations are not too fast, while the other quickly oscillates in all 
essential region of momenta y. With our choice of the sign of the initial asymmetry 
only cos A$ has an essential saddle point and in this approximation the equation ( ^T]) 
can be written as: 

y Jo 

For large q the phase difference is equal to 

A^ = K (^-^-^ + J' dq2V{q2)Z{q2)j (63) 
This integral can be taken in saddle point approximation. To this end let us expand: 
Hq, y) = Hqr, y) + ^^^^$"(gii, y) (64) 
where the saddle (resonance) point q^ is determined by the condition 

^\qR,y) = K[vZ-^ =0 (65) 

For q < qR the integral in the r.h.s. of eq. ( |0^ may be neglected, while for q> qn'ii 
is 



" dqi cos A<l> ^ eiq - qn) Tie \ J- e[*(^'?^)-*(''«-y)— /4] 1 (gg) 

Iv {"^"[QR^y)] J 

where $" = (VZ)'. 

Now repeating similar integration in eq. (^) we can easily find A(g, y): 
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Note the factor 1/2 that comes from the theta-function in expression (^). It permits 
integration only over positive values of (g — gj?). 
From the saddle point condition follows 



In the last equality the solution Z 1.5 ■ 10'^ g^^^^ and V = bq~'^^^ were used. From 
the condition {VZ)' = 1/y we find 

QR ~ (52/)' (69) 

where we used = —1 and be = 3.31 • 10^'^. 

As a simple check we may calculate the integrated lepton asymmetry L{q) = 
16.45- 10-i°Z(g): 

^, , 16.45 /•°° , 2. / \ 148i^'(sin26')2 /-y— , ,^ .q r , ^ /^,^^ 

= TFT / dyy^fe,{y)/\{q,y) = ^ / dyy {hyffe,{y) 70 

167r^ JO 167r Jo 

where ymax = This integral can be easily calculated and the result is in a good 

agreement with eq. (^9]), as one should expect. 

However the magnitude of A{q,y) is too large. For example at g = 6630 (corre- 
sponding to T = 1 MeV for {ue — z/s)-oscillations) we find A = 97c K {sin 26)'^ qj^/y ^ 
200 ^ 1. It violates the condition (^) and contradicts the assumption that A ^ S. 
Naively one would expect that the asymmetry should be suppressed by two orders of 
magnitude but as we will see below it is not the case. The evolution of Z{q) changes, 
due to the back reaction, and the behavior of the resonance yR{q) also becomes dif- 
ferent from yR = g^/^/5 found above. 

If only one resonance is essential then (S + A) is conserved, as is seen from 
eqs. (il|,il) with cosA$ —>■ 0. It corresponds to conservation of the total leptonic 
charge if oscillations are efficient only in neutrino (or antineutrino) channel. In this 
case we come to the equation: 

A'(g, y) = T dq, cos A<|. [1 - A{q,,y)] (71) 

JO 
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This integral can be taken in the same way as above in the saddle point approximation 
and we obtain: 



A(9,l/) 



qR)\[l-A{qn,y)] 



(72) 



where the last term describes the back-reaction and 

_ 7rir(sin2g)^ 
^ - f\{VZ)'\, ^''^ 

The derivative of VZ is taken over g at g = qR{y) found from the resonance condition 
V{qn)Z{qR) = 1/y. 

Since ^(0) = 1/2, we find 

2A 



A(g,i/) 



■Oiq-qn) 



2 + A 

With this expression for A{q,y) we can find integrated asymmetry 



(74) 



dyy jeq{y) 



(75) 



IQn^Jo -"^"^-^ 2 + A 

where yR^q) = l/[V{q)Z{q)]. This equation can be reduced to an ordinary differential 
equation in the following way. Let us introduce the new variable 

1 



VZ 



(76) 



and consider the new unknown function q = qij). Correspondingly Z = q'^/^ij) j {hr). 
The derivative over q should be rewritten as: 



d{VZ) 
dq 



1 /dq' 
72 [d^ 



(77) 



Under the sign of the integral over y one should take t = y, while the upper integration 
limit is ymax = t. Now we can take derivatives over r of both sides of the equation ( [75D 
and obtain: 

10^° 2A 



d 
dr 



q^l\T) 
br 



167r2 2 + A 



(78) 
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where now A = 7ri^'(sin26')^ dq/dr. This is the final equation for determination of the 
integrated asymmetry with the account of the back reaction. As initial condition we 
take the magnitude of asymmetry found from solution of eq. ( |53D at g = 5. At this q 
the back reaction is still small but already the regime of one pole dominance begins. 



Under the latter assumption the above equation (|78D is derived. 

The numerical solution of eq. ( [78l) is straightforward. In particular, if one neglects 
A with respect to 2 in the denominator of the r.h.s. of this equation, then its numerical 
solution gives exactly the same result for the asymmetry as has been found above 
from eq. (^). An account of back reaction is not essential at high temperatures 
but it is quite important in the temperature region of Big Bang Nucleosynthesis. In 
particular, at T = 1 MeV the asymmetry is L = 0.0435 and at T = 0.5 MeV it is 
L = 0.25. These values are approximately 3 times smaller than those found without 
back reaction. Asymptotic constant value of L is reached at T < 0.3 MeV and is equal 
to 0.35. These numerical values are found for electronic neutrinos with = —1 
and (sin 2^)^ = 10^^. The asymptotic value is approximately twice smaller than that 
presented in the paper the same is true for the magnitude of the asymmetry in 
the nucleosynthesis region as well. Another important effect for BBN is the shape 
of the spectrum of electronic neutrinos that may noticeably deviate from the simple 
equilibrium one given by expression @ even with a non-zero chemical potential. It 
will be considered elsewhere. 

For Un or Ur with = —10 and (sin 26')^ = 10~^ the asymmetry asymptoti- 
cally tends to 0.237 in a good agreement with ref. 0. For non-asymptotic values of 
temperature the corrected by back reaction asymmetry is = 6.94 ■ 10"'^ for T = 3 
MeV (1.24 smaller than non-corrected one), = 0.025 for T = 2 MeV (1.48 smaller), 
and Lfj_ = 0.164 for T = 1 MeV (2.6 times smaller). 

There is an easy way to find the asymptotic constant value of the asymmetry, Zq 
or Lq. To this end is convenient to use eq. (|75D with a constant Z: 

A = = 70.85(*„.V'*(Bm2«/10-)^Lr (79) 
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This result is the same both for v^, and z/^^r- 

Since the upper hmit of the integral over y is ymax = q'^^^/ (bZo) ^ 1 for large q, 
we obtain the following equation for Lq: 



Numerical solution of this equation gives L = 0.35 for Sm^ = —1 and (sin 2^)^ = 10 ^ 
and L = 0.27 for = —10 and (sin 26')^ = 10^^ in a good agreement with the 



5 Discussion. 

Thus, the statement of a generation of a large lepton asymmetry by oscillations 
between active and sterile neutrinos is essentially confirmed here. The agreement 
between the present calculations and those of ref. for the case of z/^ and mixing 
is very good, while there is a noticeable difference between the results of the present 
paper and the paper |TE[ for the case of z/g and mixing. In the temperature range 
important for the primordial nucleosynthesis the results of the present paper is 2-3 
times smaller. 

There are some other minor differences. According to our results the asymmetry 
rises as T^^/^, while the fit to numerical solution of the papers quoted in the text 
is T^. This difference is very important for the evolution of the resonance value of 
the neutrino momentum y and may possibly explain the difference of the results. On 
the other hand, if this is the case, a good agreement for the z/^ — z/^ case becomes 
mysterious. 

Another difference between the present approach and the calculations of refs. [Q- 
0,||TU|, [T^ is the treatment of repopulation of active neutrino states. In those papers 
the equilibrium number density of active neutrinos was taken in the form (^ with 
a non-zero chemical potential which value is found in a self-consistent way together 
with calculation of the lepton asymmetry. This may be true when the decoherence 




(80) 



solution of differential equation ([78|). 
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effects, described by 7, are strong. But the effects of neutrino production were impor- 
tant only at high temperatures, when a minor fraction of asymmetry was generated. 
At that stage the contribution from non-zero fi into kinetic equations are neghgible 
in comparison with the terms coming from Z in effective potential, the latter are 
amplified by the factor 10^. Moreover, if the effects of fi were non-negligible, one 
should take into account similar effects that lead to the difference between 7 and 7 
of the same magnitude. The agreement between the calculations for the (z/^ — Us) 
case discussed above shows that one may neglect fi in the equilibrium distribution 
function. 

The calculations presented here does not show any chaoticity. To be more pre- 
cise numerical solution of equation ( ^3l) shows chaotic behavior with an increasing 
i^'sin2^. However this chaoticity is related to numerical instability because with the 
increasing coefficient in front of the r.h.s. of the equation the minimal value of the 
asymmetry becomes very small and can be smaller than the accuracy of computation. 
In this case the numerically calculated value of the asymmetry may chaotically change 
sign. However this regime is well described analytically and it can be seen that the 
sign of asymmetry does not change. There could be another effect first discussed in 
ref. [|18| (see also [^) - small primordial fluctuations of the cosmological (baryonic) 
charge asymmetry could be amplified by oscillations. This effect is related to neutrino 
diffusion and has nothing to do with the discussed chaoticity. A chaotic behavior was 
observed in several recent papers [|l^, |l^, ^ in a simplified approach when the kinetic 
equations were solved for a fixed "average" value of neutrino momentum, y = 3.15, 
so that integro-differential equations become much simpler ordinary differential ones. 
However many essential features of the process could be obscured in this approach, 
in particular "running of the resonance" over neutrino spectrum and it is difficult to 
judge how reliable are the results. Moreover, the average value of 1/y that enters the 
refraction index, is (l/y) ~ 1 and not l/{y) ~ 0.3, though this numerical difference 
might not be important for the conclusion. 
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Possibly fixed-momentum approach is not adequate to the problem as can be seen 
from a very simple example. Let us consider oscillations between Ua and Ug in vacuum. 
Then for a neutrino with a fixed energy the leptonic charge in active neutrino sector 
would oscillate with a very large frequency, ~ sin{6m'^t / E) . However if one averages 
this result with thermal neutrino spectrum, the oscillations of leptonic charge would 
be exponentially suppressed. Still this counter-example is also oversimplified and 
cannot be considered as a rigorous counter- argument. It may happen for example 
that in the case of smaller K when saddle point does not give a good approximation, 
the differential asymmetry A{q,y) may be an oscillating function but the integrated 
(total) asymmetry is smooth. Another logically open option is that for a smaller K 
the asymmetry might be chaotic but not large. Momentum dependent equations were 



analyzed in ref. |jT9[ and a chaoticity for some values of parameters were observed but 
the authors could not exclude its origin by an instability of numerical computation 
procedure. At this stage is difficult to make a final judgment. 

As we have already mentioned the semi-analytic solution found in this paper 
confirms a strong generation of lepton asymmetry. The results obtained are accu- 
rate in the limit of large values of parameter K and for sufficiently small mixing, 
(sin2^^)^ < 0.01 — 0.001. For a larger mixing the population of sterile states be- 
comes very strong and the Ug and Us states closely approach equilibrium and become 
equally populated, so the asymmetry is not generated. For small values of the prod- 
uct i^' sin 26' the process of asymmetry generation is not efficient and the net result 
is rather low. Numerical calculations of the effect for very low values of the mass 
difference Sm^ = 10""^ — 10~^^ eV^ show that the asymmetry could rise only up to 4 
orders of magnitude [^-[38| producing the net result at the level of 10^^. According 
to the calculations of this work, the asymmetry strongly rises if > 10^'^ eV^ and 
possibly for smaller values depending upon the mixing angle. More accurate estimates 
of the parameter range where the asymmetry may strongly rise and the role of the 
asymmetry generation in big bang nucleosynthesis will be discussed elsewhere. 
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